Metody oparte o sąsiedztwo

Wstęp

Laboratorium składa się z 2 części:

1. kNN w uczeniu nadzorowanym
2. Wyszukiwanie

W pierwszej części użyjemy standardowych bibliotek oraz dodatkowo pynndescent, implementującego algorytm NN-Descent do approximate nearest neighbors (ANN).

W drugiej części natomiast do wektoryzacji obiektów użyjemy stosu do uczenia głębokiego opartego o PyTorcha i Sentence Transformers - biblioteki implementującej transformery do uczenia nienadzorowanego i wyszukiwania.

import matplotlib.pyplot as plt
import numpy as np
import pandas as pd

Zbiór danych - klasyfikacja danych numerycznych

Wykorzystamy zbiór danych Codon usage z dziedziny bioinformatyki. Został on zaprezentowany w artykule:

Hallee, Logan, and Bohdan B. Khomtchouk. "Machine learning classifiers predict key genomic and evolutionary traits across the kingdoms of life." Scientific Reports 13.1 (2023): 2088. link

Autorzy wykonywali na nim wiele analiz, ale nas interesuje podstawowa - przewidywanie, z jakiego królestwa (kingdom) pochodzi komórka, na podstawie rozkładu kodonów w jej sekwencjonowanym RNA. Analiza DNA oraz RNA jest podstawowym zadaniem bioinformatyki. Kodony (codons) to trójki nukleotydów, np. UGC, i większość zapisuje pewien aminokwas (z wyjątkiem trzech kodonów stopu), który jest wykorzystywany w syntezie białek.

Problem jest taki, że w praktyce trzeba sekwencjonować DNA/RNA z wielu komórek, więc trzeba je namnożyć w laboratorium. Niestety nie da się uzyskać idealnej czystości, i zawsze jest szansa, że próbka zostanie zanieczyszczona, np. bakteriami. Wtedy oprócz interesujących nas komórek (np. zwierzęcych) otrzymamy inne, których trzeba się pozbyć przed dalszą analizą. Tutaj właśnie wchodzi ML - dokonamy klasyfikacji, z jakiego królestwa (a właściwie domeny) pochodzi kod RNA, żeby pomóc w takich sytuacjach.

Dla uproszczenia, autorzy, zamiast klasyfikować królestwo, klasyfikują domenę) organizmu. Są to jednostki o poziom wyżej w systematyce taksonomicznej. Do tego dokładamy wirusy i bakteriofagi. Mamy zatem klasy:

1. Archaea (archeony) - drobne jednokomórkowce, źródło ważnych enzymów stosowanych w biologii molekularnej
2. Bacteria (bakterie)
3. Eukariota (eukarionty) - między innymi ludzie
4. Viruses (wirusy) - nie są organizmami żywymi, więc nie są częścią typowej taksonomii, ale są ważne w badaniach mikrobiologicznych
5. Bacteriophages (bakteriofagi) - rodzaj wirusów atakujących tylko bakterie, używane m.in. w badaniach nad mikrobiomem jelitowym

Szczegółowy opis zbioru znajduje się na stronie UCI. W skrócie:

1. Kingdom - królestwo, z którego pochodzi komórka. Podział jest tutaj dość szczegółowy, ale organizmy te można pogrupować w praktyce na bakterie (bacteria), wirusy (viruses) i eukarionty (eukariota).
2. DNAtype - używana w innych analizach w artykule.
3. SpeciesID - numer gatunku.
4. Ncodons - liczba kodonów zmierzona dla danej komórki.
5. SpeciesName - nazwa gatunku.

Dalsze kolumny to znormalizowany ułamek kodonów poszczególnych typów w RNA danej komórki. W związku z tym, że dane są już znormalizowane, nie ma potrzeby ich skalowania.

df = pd.read_csv("codon_usage.csv")
df

C:\Users\stas69\AppData\Local\Temp\ipykernel_10604\952778463.py:1: DtypeWarning: Columns (5,6) have mixed types. Specify dtype option on import or set low_memory=False.
  df = pd.read_csv("codon_usage.csv")

	Kingdom	DNAtype	SpeciesID	Ncodons	SpeciesName	UUU	UUC	UUA	UUG	CUU	...	CGG	AGA	AGG	GAU	GAC	GAA	GAG	UAA	UAG	UGA
0	vrl	0	100217	1995	Epizootic haematopoietic necrosis virus	0.01654	0.01203	0.00050	0.00351	0.01203	...	0.00451	0.01303	0.03559	0.01003	0.04612	0.01203	0.04361	0.00251	0.00050	0.00000
1	vrl	0	100220	1474	Bohle iridovirus	0.02714	0.01357	0.00068	0.00678	0.00407	...	0.00136	0.01696	0.03596	0.01221	0.04545	0.01560	0.04410	0.00271	0.00068	0.00000
2	vrl	0	100755	4862	Sweet potato leaf curl virus	0.01974	0.0218	0.01357	0.01543	0.00782	...	0.00596	0.01974	0.02489	0.03126	0.02036	0.02242	0.02468	0.00391	0.00000	0.00144
3	vrl	0	100880	1915	Northern cereal mosaic virus	0.01775	0.02245	0.01619	0.00992	0.01567	...	0.00366	0.01410	0.01671	0.03760	0.01932	0.03029	0.03446	0.00261	0.00157	0.00000
4	vrl	0	100887	22831	Soil-borne cereal mosaic virus	0.02816	0.01371	0.00767	0.03679	0.01380	...	0.00604	0.01494	0.01734	0.04148	0.02483	0.03359	0.03679	0.00000	0.00044	0.00131
...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...	...
13023	pri	0	9601	1097	Pongo pygmaeus abelii	0.02552	0.03555	0.00547	0.01367	0.01276	...	0.00820	0.01367	0.01094	0.01367	0.02279	0.02005	0.04102	0.00091	0.00091	0.00638
13024	pri	1	9601	2067	mitochondrion Pongo pygmaeus abelii	0.01258	0.03193	0.01984	0.00629	0.01451	...	0.00145	0.00000	0.00048	0.00194	0.01306	0.01838	0.00677	0.00242	0.00097	0.01887
13025	pri	1	9602	1686	mitochondrion Pongo pygmaeus pygmaeus	0.01423	0.03321	0.01661	0.00356	0.01127	...	0.00000	0.00000	0.00000	0.00178	0.01661	0.02788	0.00297	0.00356	0.00119	0.02017
13026	pri	0	9606	40662582	Homo sapiens	0.01757	0.02028	0.00767	0.01293	0.01319	...	0.01142	0.01217	0.01196	0.02178	0.02510	0.02896	0.03959	0.00099	0.00079	0.00156
13027	pri	1	9606	8998998	mitochondrion Homo sapiens	0.01778	0.03724	0.01732	0.00600	0.01689	...	0.00083	0.00041	0.00041	0.00451	0.01402	0.01651	0.00783	0.00156	0.00114	0.02161

13028 rows × 69 columns

Mamy warning co do typu - niedobrze. Sprawdźmy typy:

df.dtypes.head(20)

Kingdom         object
DNAtype          int64
SpeciesID        int64
Ncodons          int64
SpeciesName     object
UUU             object
UUC             object
UUA            float64
UUG            float64
CUU            float64
CUC            float64
CUA            float64
CUG            float64
AUU            float64
AUC            float64
AUA            float64
AUG            float64
GUU            float64
GUC            float64
GUA            float64
dtype: object

Kolumny UUU i UUC powinny być floatami - niedobrze. Na szczęście to tylko literówka w jednym wierszu, co można sprawdzić w danych. Usuniemy ją po prostu.

df = df[pd.to_numeric(df["UUU"], errors="coerce").notnull()].copy()

df = df.copy()  # to avoid irritating SettingWithCopyWarning

df["UUU"] = df.loc[:, "UUU"].astype(float)
df["UUC"] = df.loc[:, "UUC"].astype(float)

df.dtypes.head(20)

Kingdom         object
DNAtype          int64
SpeciesID        int64
Ncodons          int64
SpeciesName     object
UUU            float64
UUC            float64
UUA            float64
UUG            float64
CUU            float64
CUC            float64
CUA            float64
CUG            float64
AUU            float64
AUC            float64
AUA            float64
AUG            float64
GUU            float64
GUC            float64
GUA            float64
dtype: object

Zadanie 1 (1.5 punktu)

1. Usuń wiersze mające mniej niż 1000 kodonów
2. Usuń wiersze dla plazmidów (wartość plm dla królestwa)
3. Usuń kolumny bezużyteczne dla ML: DNAtype, SpeciesID, Ncodons, SpeciesName
4. Zakoduj klasy (kolumna Kingdom):

• 0 - archaea, arc
• 1 - bacteria, bct
• 2 - eukaryota, pln, inv, vrt, mam, rod, pri
• 3 - bacteriophages (phages), phg
• 4 - viruses (viral cells), vrl

5. Wyodrębnij klasy jako osobną zmienną y. Pamiętaj o usunięciu tej kolumny z oryginalnego DataFrame'a.
6. Narysuj wykres rozkładu klas:
   • wykres słupkowy (bar plot)
   • na osi X mają być nazwy klas pod słupkami, wypisane poziomo
   • wartości na osi Y mają być w procentach (znormalizowane 0-100)
   • pamiętaj o tytule oraz opisaniu osi

# 1
df = df[df["Ncodons"] >= 1000]
# 2
df = df[df["Kingdom"] != "plm"]
# 3
df.drop(columns=["DNAtype", "SpeciesID", "Ncodons", "SpeciesName"], inplace=True)
# 4
kingdom_mapping = {
    "arc": 0,
    "bct": 1,
    "pln": 2,
    "inv": 2,
    "vrt": 2,
    "mam": 2,
    "rod": 2,
    "pri": 2,
    "phg": 3,
    "vrl": 4,
}
df["Kingdom"] = df["Kingdom"].map(kingdom_mapping)

# 5
y = df["Kingdom"]
df.drop(columns=["Kingdom"], inplace=True)

# 6
kigdom_mapping = {
    0: "archaea",
    1: "bacteria",
    2: "eukaryota",
    3: "bacteriophages",
    4: "viruses",
}
class_counts = y.value_counts(normalize=True) * 100
class_counts.index = class_counts.index.map(kigdom_mapping)
bar = class_counts.plot.bar()
for i, v in enumerate(class_counts):
    bar.text(i, v + 1, f"{v:.1f}%", ha="center", fontsize=10)

plt.xticks(rotation=0)
plt.xlabel("Class (Kingdom)")
plt.ylabel("Percentage [%]")
plt.title("Class Percentage in Dataset")

Text(0.5, 1.0, 'Class Percentage in Dataset')

Jak widać, mamy nie tylko klasyfikację wieloklasową, ale na dodatek klasyfikację niezbalansowaną. W takim przypadku trzeba użyć odpowiednich metryk. Wzorując się na artykule, użyjemy macro-averaged F1-score. Jest to proste rozwinięcie F1-score na wiele klas.

Dla przypomnienia, w klasyfikacji binarnej, F1-score to średnia harmoniczna precyzji (precision) i czułości (recall):

$$ F1 = \frac{2}{recall^{-1} + precision^{-1}} = 2 * \frac{precision * recall}{precision + recall} = \frac{2 * TP}{2 * TP + FP + FN} $$

Macro-averaging oznacza, że najpierw obliczamy F1-score dla każdej klasy, a później bierzemy ich średnią. Żeby taka średnia była wysoka, to musimy mieć wysoki wynik na wszystkich 5 klasach, w tym 2 mniejszościowych, więc uwzględnimy niezbalansowanie klas.

Klasyfikator kNN

Wytrenujemy teraz podstawowy algorytm k najbliższych sąsiadów, korzystając z samych wartości domyślnych ze Scikit-learn: 5 sąsiadów, metryka euklidesowa, brak ważenia.

Zgodnie z artykułem, wykorzystamy zwykłą metodę holdout do testowania, dzieląc zbiór na treningowy i testowy w proporcjach 80%-20%. Zrobimy to ze stratyfikacją (stratification), czyli tak, żeby rozkład klas był taki sam w obu podzbiorach.

from sklearn.model_selection import train_test_split


X_train, X_test, y_train, y_test = train_test_split(
    df, y, test_size=0.2, random_state=0, stratify=y
)

from sklearn.metrics import f1_score
from sklearn.neighbors import KNeighborsClassifier


clf = KNeighborsClassifier(n_jobs=-1)
clf.fit(X_train, y_train)

y_pred_score_train = clf.predict(X_train)
y_pred_score_test = clf.predict(X_test)

f1_train = f1_score(y_train, y_pred_score_train, average="macro")
f1_test = f1_score(y_test, y_pred_score_test, average="macro")

print(f"F1 train: {100 * f1_train:.2f}%")
print(f"F1 test: {100 * f1_test:.2f}%")

F1 train: 92.51%
F1 test: 91.35%

Wyniki są całkiem niezłe, ale jest jeszcze miejsce na poprawę. Być może zmiana liczby sąsiadów albo metryki, albo też dodanie ważenia sąsiadów da nam jeszcze kilka procent.

W przypadku wielu klas samo agregowane F1 nie zawsze mówi jednak wszystko. Zawsze warto sprawdzić macierz pomyłek (confusion matrix). Co ważne, można ją rysować albo z surową liczbą przykładów, albo znormalizować, żeby dostać wyniki w ułamku (zakres wartości 0-1). Pierwsza możliwość jest przydatna, kiedy mamy niezbalansowane klasy i patrzymy na ogólną jakość klasyfikatora. Znormalizowane wykresy są z kolei przydatne, kiedy przyglądamy się klasom z osobna.

from sklearn.metrics import ConfusionMatrixDisplay

ConfusionMatrixDisplay.from_estimator(clf, X_test, y_test)
plt.title("Confusion matrix (no normalization)")
plt.show()


ConfusionMatrixDisplay.from_estimator(clf, X_test, y_test, normalize="true")
plt.title("Confusion matrix (row-normalized)")
plt.show()

Wyniki są całkiem niezłe, w końcu mamy wysokie F1. Widać jednak, że dość często mylimy klasę 0 (archeony) oraz klasę 3 (bakteriofagi) z klasą 1 (bakterie) - jest to problematyczne, ale zrozumiałe, bo tych klas jest najmniej.

Po tuningu powinniśmy otrzymać lepsze wyniki. Na dobry początek sprawdzimy różną liczbę sąsiadów. Trzeba zwrócić uwagę, że jeżeli nasza metryka ma jakieś argumenty, to trzeba użyć funkcji make_scorer i przekazać wynik do GridSearchCV.

Zwróć też uwagę, że skoro przekazujemy n_jobs=-1 do klasyfikatora, to grid search dostaje n_jobs=None, żeby mieć tyle procesów, co rdzeni procesora. Przy okazji zmierzymy też czas tuningu.

from time import time

from sklearn.metrics import make_scorer
from sklearn.model_selection import GridSearchCV


clf = KNeighborsClassifier(n_jobs=-1)

param_grid = {
    "n_neighbors": list(range(1, 51)),
}

multiclass_f1 = make_scorer(
    f1_score,
    average="macro",
    greater_is_better=True,
)

cv = GridSearchCV(
    estimator=clf,
    param_grid=param_grid,
    scoring=multiclass_f1,
    cv=5,
)

time_start = time()
cv.fit(X_train, y_train)
time_end = time()

print(f"Optimal number of neighbors: {cv.best_params_['n_neighbors']}")
print(f"Tuning time: {time_end - time_start:.2f} s")

y_pred_score_train = cv.predict(X_train)
y_pred_score_test = cv.predict(X_test)

f1_train = f1_score(y_train, y_pred_score_train, average="macro")
f1_test = f1_score(y_test, y_pred_score_test, average="macro")

print(f"F1 train: {100 * f1_train:.2f}%")
print(f"F1 test: {100 * f1_test:.2f}%")

Optimal number of neighbors: 1
Tuning time: 31.11 s
F1 train: 100.00%
F1 test: 89.16%

Co ciekawe, po tuningu nie tylko nie mamy lepszego wyniku, ale wręcz przeuczamy bardziej! Tak się niestety zdarza - tuning opiera się na wyniku na zbiorze walidacyjnym, które nie zawsze są dobrą estymatą. Pewnym sposobem na ominięcie tego jest zmniejszenie zakresu dla hiperparametru. Tutaj wiemy, że:

• im mniejsza liczba sąsiadów, tym większa wariancja (mocniejszy overfitting) i na odwrót
• domyślna liczba sąsiadów to 5
• przy 4 przeuczamy mocniej

Można więc po prostu sprawdzać zakres od 5 wzwyż. Do tego czeka nas jeszcze tuning metryki i ważenia sąsiadów.

Zadanie 2 (1.5 punktu)

1. Dokonaj tuningu: liczby sąsiadów (zakres [5, 25]), metryki (euklidesowa, Manhattan, cosinusowa) oraz ważenia sąsiadów (brak vs ważenie odległością). Jakie są optymalne wartości hiperparametrów?
2. Zmierz czas. Czy twoim zdaniem to długo, dla zbioru, który ma ok. 10k próbek uczących?
3. Narysuj macierz pomyłek dla zbioru testowgeo. Czy wygląda sensownie, szczególnie dla klas mniejszościowych? Czy udało się uzyskać wynik lepszy od bazowego?

clf = KNeighborsClassifier(n_jobs=-1)

param_grid = {
    "n_neighbors": list(range(5, 26)),
    "metric": ["euclidean", "cosine", "manhattan"],
    "weights": ["uniform", "distance"],
}

multiclass_f1 = make_scorer(
    f1_score,
    average="macro",
    greater_is_better=True,
)

cv = GridSearchCV(
    estimator=clf,
    param_grid=param_grid,
    scoring=multiclass_f1,
    cv=5,
)

time_start = time()
cv.fit(X_train, y_train)
time_end = time()

print(f"Optimal number of neighbors: {cv.best_params_['n_neighbors']}")
print(f"Tuning time: {time_end - time_start:.2f} s")

y_pred_score_train = cv.predict(X_train)
y_pred_score_test = cv.predict(X_test)

f1_train = f1_score(y_train, y_pred_score_train, average="macro")
f1_test = f1_score(y_test, y_pred_score_test, average="macro")

print(f"F1 train: {100 * f1_train:.2f}%")
print(f"F1 test: {100 * f1_test:.2f}%")

Optimal number of neighbors: 6
Tuning time: 136.35 s
F1 train: 100.00%
F1 test: 91.53%

print("Best Metric:", cv.best_params_["metric"])
print("Best Weights:", cv.best_params_["weights"])

Best Metric: cosine
Best Weights: distance

Wydaje mi się, że 3 minuty to nie tak długo szczególnie patrząc na liczbę próbek, a także na grid searcha z poprzedniego labu.

from sklearn.metrics import ConfusionMatrixDisplay

best_clf = cv.best_estimator_

best_clf.fit(X_train, y_train)
ConfusionMatrixDisplay.from_estimator(best_clf, X_test, y_test)
plt.title("Confusion matrix (no normalization)")
plt.show()


ConfusionMatrixDisplay.from_estimator(best_clf, X_test, y_test, normalize="true")
plt.title("Confusion matrix (row-normalized)")
plt.show()

Udało się przebić model bazowy o 0,18%. Macierz jest w zasadzie bardzo podobna do tej poprzedniej.

Approximate Nearest Neighbors (ANN)

Postaramy się teraz przyspieszyć nasz klasyfikator, w zamian za odrobinę dokładności, za pomocą algorytmu przybliżonych najbliższych sąsiadów, a konkretnie NN-Descent. Ten został wybrany głównie ze względu na łatwą w użyciu implementację.

Twórcy Scikit-learn'a zdawali sobie sprawę, że algorytmy znajdowania najbliższych sąsiadów to aktywne pole rozwoju, dlatego stworzyli klasę (i interfejs) KNeighborsTransformer, pozwalającą integrować zewnętrzne rozwiązania. Działa ona zarówno do klasyfikacji / regresji kNN, jak i do innych algorytmów opartych o sąsiedztwo, które budują graf najbliższych sąsiadów. W przypadku jego użycia trzeba podać metric="precomputed", bo obliczaniem odległości zajmuje się osobny algorytm.

Implementacja PyNNDescent korzysta właśnie z tego rozwiązania, co gwarantuje łatwość użycia - klasa PyNNDescentTransformer dziedziczy po KNeighborsTransformer. Zobaczmy, jak to działa. Wykorzystamy funkcję make_pipeline(), która buduje obiekt Pipeline, ale nie wymaga podawania nazw kolejnych elementów.

Uwaga: PyNNDescent dla wydajności korzysta z Numby, czyli kompilatora JIT (Just-In-Time). Dzięki temu działa to bardzo szybko, ale pierwszy import i wykonanie zajmuje sporo czasu. Wykonaj kod 2 razy - za drugim razem będzie o wiele szybciej.

from time import time


def benchmark_knn_and_ann(
    sklearn_knn,
    pynndescent_ann,
    X_train,
    X_test,
    y_train,
    y_test,
) -> None:
    # training
    start_time = time()
    sklearn_knn.fit(X_train, y_train)
    end_time = time()

    sklearn_knn_fit_time = end_time - start_time

    start_time = time()
    pynndescent_ann.fit(X_train, y_train)
    end_time = time()

    pynndescent_ann_fit_time = end_time - start_time

    # prediction

    start_time = time()
    y_pred_sklearn = sklearn_knn.predict(X_test)
    end_time = time()

    sklearn_knn_predict_time = end_time - start_time

    start_time = time()
    y_pred_pynndescent = pynndescent_ann.predict(X_test)
    end_time = time()

    pynndescent_ann_predict_time = end_time - start_time

    f1_knn = f1_score(y_test, y_pred_sklearn, average="macro")
    f1_ann = f1_score(y_test, y_pred_pynndescent, average="macro")

    print(f"Scikit-learn training time: {sklearn_knn_fit_time:.2f}")
    print(f"PyNNDescent training time: {pynndescent_ann_fit_time:.2f}")
    print()
    print(f"Scikit-learn prediction time: {sklearn_knn_predict_time:.2f}")
    print(f"PyNNDescent prediction time: {pynndescent_ann_predict_time:.2f}")
    print()
    print(f"Scikit-learn F1: {100 * f1_knn:.2f}%")
    print(f"PyNNDescent F1: {100 * f1_ann:.2f}%")

from pynndescent import PyNNDescentTransformer
from sklearn.pipeline import make_pipeline


sklearn_knn = KNeighborsClassifier(metric="euclidean")

pynndescent_ann = make_pipeline(
    PyNNDescentTransformer(metric="euclidean", random_state=0),
    KNeighborsClassifier(metric="precomputed"),
)

benchmark_knn_and_ann(sklearn_knn, pynndescent_ann, X_train, X_test, y_train, y_test)

Scikit-learn training time: 0.00
PyNNDescent training time: 25.95

Scikit-learn prediction time: 0.14
PyNNDescent prediction time: 0.14

Scikit-learn F1: 91.35%
PyNNDescent F1: 91.25%

Czas treningu jest dłuższy - to oczekiwany wynik. Natomiast czas predykcji jest podobny - wynika to z tego, że dla metryki euklidesowej k-d tree jest bardzo wydajne. Wynik jest bardzo podobny, więc mamy dobre przybliżenie.

Sprawdźmy metrykę Manhattan, dla której typowo k-d tree radzi sobie nieco gorzej.

sklearn_knn = KNeighborsClassifier(metric="manhattan")

pynndescent_ann = make_pipeline(
    PyNNDescentTransformer(metric="manhattan", random_state=0),
    KNeighborsClassifier(metric="precomputed"),
)

benchmark_knn_and_ann(sklearn_knn, pynndescent_ann, X_train, X_test, y_train, y_test)

Scikit-learn training time: 0.01
PyNNDescent training time: 10.77

Scikit-learn prediction time: 0.28
PyNNDescent prediction time: 0.15

Scikit-learn F1: 91.39%
PyNNDescent F1: 91.39%

Tutaj predykcja ANN jest już wyraźnie szybsza - co prawda dane są niezbyt duże, więc nie ma to aż takiego znaczenia, ale widać potencjał skalowalności. Dodatkowo wynik jest taki sam.

Zbiór danych - klasyfikacja danych chemicznych

Teraz zajmiemy się klasyfikacją zbioru posiadającego same zmienne binarne, do których trzeba użyć odpowiednich metryk. Konkretnie, będzie to zbiór HIV z benchmarku MoleculeNet, w którym na podstawie reprezentacji chemicznej molekuły trzeba przewidywać, czy jest ona inhibitorem wirusa HIV. Inhibitory to substancje przeciwne do katalizatorów - spowalniają lub uniemożliwiają jakieś zjawisko, w naszym wypadku infekcję wirusa HIV. Przewidywanie własności molekuł to kluczowe zadanie w nowoczesnym projektowaniu leków i bardzo ważne zastosowanie ML w farmacji i chemii molekularnej.

Zbiór jest już podzielony na treningowy, walidacyjny i testowy w ramach projektu Open Graph Benchmark. Dane te można traktować jak klasyfikację grafów molekularnych, np. z pomocą grafowych sieci neuronowych. Jednak podejściem klasycznym i dającym często lepsze wyniki są fingerprinty molekularne (molecular fingerprints). Są to deterministyczne algorytmy ekstrakcji cech z grafu, zamieniające go na wektor. Jest ich dużo, ale najpopularniejszy to ECFP (Extended Connectivity FingerPrint), wykorzystujący informację o subgrafach o małym promieniu (typowo 4). Domyślnie skutkuje to cechami binarnymi, czy subgraf wystąpił, czy nie w danej molekule.

Molekuły są typowo przechowywane w formacie SMILES strings. Graf molekularny można odczytać z tego formatu, a potem przekazać taki zbiór do fingerprintu. Wtedy zamieniamy problem na zwykłą klasyfikację tabelaryczną. Scikit-fingerprints zrobi to automatycznie, jeżeli SMILES strings zostaną przekazane do fingerprintu.

Jako że zbiór jest binarny, to nie ma potrzeby żadnego skalowania zmiennych. W przypadku zbioru HIV typową metryką jest Area Under Receiver Operating Characteristic (AUROC / ROC AUC).

Zadanie 3 (1 punkt)

1. Załaduj zbiór HIV z pomocą funkcji load_hiv(). Wypisz kilka pierwszych SMILESów.
2. Załaduj indeksy treningowe, walidacyjne i testowe z benchmarku OGB z pomocą funkcji load_ogb_splits().
3. Podziel zbiór (SMILES oraz klasy y) na treningowe, walidacyjne i treningowe.
4. Narysuj wykres częstości klas w zbiorze treningowym. Jak sądzisz, czemu akurat AUROC zostało tutaj wybrane jako metryka?
5. Oblicz fingerprinty ECFP z użyciem ECFPFingerprint dla wszystkich 3 podzbiorów.

from skfp.datasets.moleculenet import load_hiv, load_ogb_splits

dataset = load_hiv(as_frame=True)
dataset.head()

C:\Users\stas69\Documents\Studia\podstawy-uczenia-maszynowego-24-25\.venv\Lib\site-packages\huggingface_hub\file_download.py:797: FutureWarning: `resume_download` is deprecated and will be removed in version 1.0.0. Downloads always resume when possible. If you want to force a new download, use `force_download=True`.
  warnings.warn(

	SMILES	label
0	CCC1=[O+][Cu-3]2([O+]=C(CC)C1)[O+]=C(CC)CC(CC)...	0
1	C(=Cc1ccccc1)C1=[O+][Cu-3]2([O+]=C(C=Cc3ccccc3...	0
2	CC(=O)N1c2ccccc2Sc2c1ccc1ccccc21	0
3	Nc1ccc(C=Cc2ccc(N)cc2S(=O)(=O)O)c(S(=O)(=O)O)c1	0
4	O=S(=O)(O)CCS(=O)(=O)O	0

split_indexes = load_ogb_splits(dataset_name="HIV", as_dict=True)
X = dataset["SMILES"]
y = dataset["label"]
X_train = X.iloc[split_indexes["train"]]
y_train = y.iloc[split_indexes["train"]]

X_test = X.iloc[split_indexes["test"]]
y_test = y.iloc[split_indexes["test"]]

X_val = X.iloc[split_indexes["valid"]]
y_val = y.iloc[split_indexes["valid"]]

C:\Users\stas69\Documents\Studia\podstawy-uczenia-maszynowego-24-25\.venv\Lib\site-packages\huggingface_hub\file_download.py:797: FutureWarning: `resume_download` is deprecated and will be removed in version 1.0.0. Downloads always resume when possible. If you want to force a new download, use `force_download=True`.
  warnings.warn(

class_counts = y_train.value_counts(normalize=True) * 100
bar = class_counts.plot.bar()
for i, v in enumerate(class_counts):
    bar.text(i, v + 1, f"{v:.1f}%", ha="center", fontsize=10)

plt.xticks(rotation=0)
plt.xlabel("Classes")
plt.ylabel("Percentage [%]")
plt.title("Class Percentage in y_train")

Text(0.5, 1.0, 'Class Percentage in y_train')

Jak widać na wykresie, zbiór danych jest mocno niezbalansowany, klasa 0 występuje o wiele częściej niż 1. W takiej sytuacji metryka AUROC jest trafnym wyborem do oceny jakości modelu. W przeciwieństwie do accuracy, która po prostu mierzy odsetek poprawnych klasyfikacji (i może być sztucznie zawyżona przez przewidywanie klasy większościowej), AUROC ocenia zdolność modelu do rozróżniania pomiędzy klasami niezależnie od ich proporcji. Dzięki temu jest bardziej odporna na problemy wynikające z niezbalansowanego rozkładu klas i lepiej oddaje rzeczywistą skuteczność klasyfikatora.

from skfp.fingerprints import ECFPFingerprint

ecfp = ECFPFingerprint()
X_train_fp = ecfp.transform(X_train)
X_val_fp = ecfp.transform(X_val)
X_test_fp = ecfp.transform(X_test)

[10:27:31] WARNING: not removing hydrogen atom without neighbors
[10:27:31] WARNING: not removing hydrogen atom without neighbors

Zadanie 4 (1.5 punktu)

1. Wytrenuj klasyfikator kNN z domyślnymi hiperparametrami, używając metryki Jaccarda (Tanimoto). Zmierz AUROC na zbiorze testowym. Pamiętaj, żeby przekazać odpowiednie wartości do metryki (patrz przykład w dokumentacji).
2. Zmierz czas treningu oraz predykcji. Czy podczas treningu została zbudowana jakakolwiek struktura danych?
3. Dokonaj tuningu liczby sąsiadów z zakresu [1, 5, 10]. Wykorzystaj tutaj istniejący zbiór walidacyjny z pomocą PredefinedSplit. Może się przydać ta odpowiedź na StackOverflow. Do tuningu wykorzystaj metrykę AUROC (argument scoring).
4. Porównaj wyniki algorytmu po tuningu z wynikami z artykułu "Molecular Topological Profile (MOLTOP) -- Simple and Strong Baseline for Molecular Graph Classification" J. Adamczyk, W. Czech.

clf = KNeighborsClassifier(n_jobs=-1, metric="jaccard")

%%time
clf.fit(X_train_fp, y_train)

CPU times: total: 31.2 ms
Wall time: 42 ms

KNeighborsClassifier(metric='jaccard', n_jobs=-1)

In a Jupyter environment, please rerun this cell to show the HTML representation or trust the notebook.
On GitHub, the HTML representation is unable to render, please try loading this page with nbviewer.org.

%%time
y_test_predict = clf.predict(X_test_fp)

C:\Users\stas69\Documents\Studia\podstawy-uczenia-maszynowego-24-25\.venv\Lib\site-packages\sklearn\metrics\pairwise.py:2466: DataConversionWarning: Data was converted to boolean for metric jaccard
  warnings.warn(msg, DataConversionWarning)
C:\Users\stas69\Documents\Studia\podstawy-uczenia-maszynowego-24-25\.venv\Lib\site-packages\sklearn\metrics\pairwise.py:2466: DataConversionWarning: Data was converted to boolean for metric jaccard
  warnings.warn(msg, DataConversionWarning)

CPU times: total: 14min 39s
Wall time: 1min 16s

from sklearn.metrics import roc_auc_score

auroc = roc_auc_score(y_test, y_test_predict, average="macro")
print(f"AUROC: {100 * auroc:.2f}%")

AUROC: 59.02%

from sklearn.model_selection import PredefinedSplit
from sklearn.model_selection import GridSearchCV

param_grid = {"n_neighbors": [1, 5, 10]}

X_train_val = np.concatenate((X_train_fp, X_val_fp), axis=0)
y_train_val = np.concatenate((y_train, y_val), axis=0)

split_index = [-1] * len(X_train_fp) + [0] * len(X_val_fp)
ps = PredefinedSplit(test_fold=split_index)

knn = KNeighborsClassifier(n_jobs=-1, metric="jaccard")

grid = GridSearchCV(estimator=knn, param_grid=param_grid, cv=ps, scoring="roc_auc")
grid.fit(X_train_val, y_train_val)


print("Best params:", grid.best_params_)
print(f"Best AUROC (val): {100 * grid.best_score_:.2f}%")

C:\Users\stas69\Documents\Studia\podstawy-uczenia-maszynowego-24-25\.venv\Lib\site-packages\sklearn\metrics\pairwise.py:2466: DataConversionWarning: Data was converted to boolean for metric jaccard
  warnings.warn(msg, DataConversionWarning)
C:\Users\stas69\Documents\Studia\podstawy-uczenia-maszynowego-24-25\.venv\Lib\site-packages\sklearn\metrics\pairwise.py:2466: DataConversionWarning: Data was converted to boolean for metric jaccard
  warnings.warn(msg, DataConversionWarning)
C:\Users\stas69\Documents\Studia\podstawy-uczenia-maszynowego-24-25\.venv\Lib\site-packages\sklearn\metrics\pairwise.py:2466: DataConversionWarning: Data was converted to boolean for metric jaccard
  warnings.warn(msg, DataConversionWarning)
C:\Users\stas69\Documents\Studia\podstawy-uczenia-maszynowego-24-25\.venv\Lib\site-packages\sklearn\metrics\pairwise.py:2466: DataConversionWarning: Data was converted to boolean for metric jaccard
  warnings.warn(msg, DataConversionWarning)
C:\Users\stas69\Documents\Studia\podstawy-uczenia-maszynowego-24-25\.venv\Lib\site-packages\sklearn\metrics\pairwise.py:2466: DataConversionWarning: Data was converted to boolean for metric jaccard
  warnings.warn(msg, DataConversionWarning)
C:\Users\stas69\Documents\Studia\podstawy-uczenia-maszynowego-24-25\.venv\Lib\site-packages\sklearn\metrics\pairwise.py:2466: DataConversionWarning: Data was converted to boolean for metric jaccard
  warnings.warn(msg, DataConversionWarning)

Best params: {'n_neighbors': 10}
Best AUROC (val): 74.47%

y_pred_cv_train = grid.predict(X_train_fp)
y_pred_cv_test = grid.predict(X_test_fp)
y_pred_cv_val = grid.predict(X_val_fp)

print(f"Best AUROC (train): {100 * roc_auc_score(y_train, y_pred_cv_train):.2f}%")
print(f"Best AUROC (test): {100 * roc_auc_score(y_test, y_pred_cv_test):.2f}%")
print(f"Best AUROC (val): {100 * roc_auc_score(y_val, y_pred_cv_val):.2f}%")

C:\Users\stas69\Documents\Studia\podstawy-uczenia-maszynowego-24-25\.venv\Lib\site-packages\sklearn\metrics\pairwise.py:2466: DataConversionWarning: Data was converted to boolean for metric jaccard
  warnings.warn(msg, DataConversionWarning)
C:\Users\stas69\Documents\Studia\podstawy-uczenia-maszynowego-24-25\.venv\Lib\site-packages\sklearn\metrics\pairwise.py:2466: DataConversionWarning: Data was converted to boolean for metric jaccard
  warnings.warn(msg, DataConversionWarning)
C:\Users\stas69\Documents\Studia\podstawy-uczenia-maszynowego-24-25\.venv\Lib\site-packages\sklearn\metrics\pairwise.py:2466: DataConversionWarning: Data was converted to boolean for metric jaccard
  warnings.warn(msg, DataConversionWarning)
C:\Users\stas69\Documents\Studia\podstawy-uczenia-maszynowego-24-25\.venv\Lib\site-packages\sklearn\metrics\pairwise.py:2466: DataConversionWarning: Data was converted to boolean for metric jaccard
  warnings.warn(msg, DataConversionWarning)
C:\Users\stas69\Documents\Studia\podstawy-uczenia-maszynowego-24-25\.venv\Lib\site-packages\sklearn\metrics\pairwise.py:2466: DataConversionWarning: Data was converted to boolean for metric jaccard
  warnings.warn(msg, DataConversionWarning)
C:\Users\stas69\Documents\Studia\podstawy-uczenia-maszynowego-24-25\.venv\Lib\site-packages\sklearn\metrics\pairwise.py:2466: DataConversionWarning: Data was converted to boolean for metric jaccard
  warnings.warn(msg, DataConversionWarning)
C:\Users\stas69\Documents\Studia\podstawy-uczenia-maszynowego-24-25\.venv\Lib\site-packages\sklearn\metrics\pairwise.py:2466: DataConversionWarning: Data was converted to boolean for metric jaccard
  warnings.warn(msg, DataConversionWarning)
C:\Users\stas69\Documents\Studia\podstawy-uczenia-maszynowego-24-25\.venv\Lib\site-packages\sklearn\metrics\pairwise.py:2466: DataConversionWarning: Data was converted to boolean for metric jaccard
  warnings.warn(msg, DataConversionWarning)
C:\Users\stas69\Documents\Studia\podstawy-uczenia-maszynowego-24-25\.venv\Lib\site-packages\sklearn\metrics\pairwise.py:2466: DataConversionWarning: Data was converted to boolean for metric jaccard
  warnings.warn(msg, DataConversionWarning)
C:\Users\stas69\Documents\Studia\podstawy-uczenia-maszynowego-24-25\.venv\Lib\site-packages\sklearn\metrics\pairwise.py:2466: DataConversionWarning: Data was converted to boolean for metric jaccard
  warnings.warn(msg, DataConversionWarning)
C:\Users\stas69\Documents\Studia\podstawy-uczenia-maszynowego-24-25\.venv\Lib\site-packages\sklearn\metrics\pairwise.py:2466: DataConversionWarning: Data was converted to boolean for metric jaccard
  warnings.warn(msg, DataConversionWarning)
C:\Users\stas69\Documents\Studia\podstawy-uczenia-maszynowego-24-25\.venv\Lib\site-packages\sklearn\metrics\pairwise.py:2466: DataConversionWarning: Data was converted to boolean for metric jaccard
  warnings.warn(msg, DataConversionWarning)
C:\Users\stas69\Documents\Studia\podstawy-uczenia-maszynowego-24-25\.venv\Lib\site-packages\sklearn\metrics\pairwise.py:2466: DataConversionWarning: Data was converted to boolean for metric jaccard
  warnings.warn(msg, DataConversionWarning)
C:\Users\stas69\Documents\Studia\podstawy-uczenia-maszynowego-24-25\.venv\Lib\site-packages\sklearn\metrics\pairwise.py:2466: DataConversionWarning: Data was converted to boolean for metric jaccard
  warnings.warn(msg, DataConversionWarning)

Best AUROC (train): 69.53%
Best AUROC (test): 57.96%
Best AUROC (val): 66.54%

Wyniki otrzymane tutaj są o wiele niższe niż te przedstawione w artykule tam wynoszą około 80% dla zbioru hiv. Wydaje mi się że nie została zbudowana żadna struktura co może świadczyć o długim czasie predykcji.

Zadanie 5 (1.5 punktu)

1. Stwórz pipeline do klasyfikacji ANN z użyciem PyNNDescent. Wykorzystaj metrykę Jaccarda (Tanimoto) i znalezioną wcześniej optymalną liczbę sąsiadów.
2. Wytrenuj zwykłe kNN oraz ANN oraz dokonaj predykcji. Zmierz czas oraz dokładność obu metod. Czy twoim zdaniem warto dokonać w tym wypadku takiej aproksymacji?

pynndescent_ann = make_pipeline(
    PyNNDescentTransformer(n_neighbors=10, metric="jaccard", random_state=0),
    KNeighborsClassifier(metric="precomputed", n_neighbors=10),
)
sklearn_knn = KNeighborsClassifier(metric="jaccard", n_neighbors=10, n_jobs=-1)

benchmark_knn_and_ann(
    sklearn_knn, pynndescent_ann, X_train_fp, X_test_fp, y_train, y_test
)

C:\Users\stas69\Documents\Studia\podstawy-uczenia-maszynowego-24-25\.venv\Lib\site-packages\sklearn\metrics\pairwise.py:2466: DataConversionWarning: Data was converted to boolean for metric jaccard
  warnings.warn(msg, DataConversionWarning)
C:\Users\stas69\Documents\Studia\podstawy-uczenia-maszynowego-24-25\.venv\Lib\site-packages\sklearn\metrics\pairwise.py:2466: DataConversionWarning: Data was converted to boolean for metric jaccard
  warnings.warn(msg, DataConversionWarning)

Scikit-learn training time: 0.00
PyNNDescent training time: 5.66

Scikit-learn prediction time: 77.72
PyNNDescent prediction time: 0.68

Scikit-learn F1: 61.36%
PyNNDescent F1: 60.72%

from sklearn.metrics import roc_auc_score


def benchmark_knn_and_ann_auc(
    sklearn_knn,
    pynndescent_ann,
    X_train,
    X_test,
    y_train,
    y_test,
) -> None:
    # training
    start_time = time()
    sklearn_knn.fit(X_train, y_train)
    end_time = time()

    sklearn_knn_fit_time = end_time - start_time

    start_time = time()
    pynndescent_ann.fit(X_train, y_train)
    end_time = time()

    pynndescent_ann_fit_time = end_time - start_time

    # prediction

    start_time = time()
    y_pred_sklearn = sklearn_knn.predict(X_test)
    end_time = time()

    sklearn_knn_predict_time = end_time - start_time

    start_time = time()
    y_pred_pynndescent = pynndescent_ann.predict(X_test)
    end_time = time()

    pynndescent_ann_predict_time = end_time - start_time

    roc_auc_knn = roc_auc_score(y_test, y_pred_sklearn, average="macro")
    roc_auc_ann = roc_auc_score(y_test, y_pred_pynndescent, average="macro")

    print(f"Scikit-learn training time: {sklearn_knn_fit_time:.2f}")
    print(f"PyNNDescent training time: {pynndescent_ann_fit_time:.2f}")
    print()
    print(f"Scikit-learn prediction time: {sklearn_knn_predict_time:.2f}")
    print(f"PyNNDescent prediction time: {pynndescent_ann_predict_time:.2f}")
    print()
    print(f"Scikit-learn AUROC: {100 * roc_auc_knn:.2f}%")
    print(f"PyNNDescent AUROC: {100 * roc_auc_ann:.2f}%")

benchmark_knn_and_ann_auc(
    sklearn_knn, pynndescent_ann, X_train_fp, X_test_fp, y_train, y_test
)

C:\Users\stas69\Documents\Studia\podstawy-uczenia-maszynowego-24-25\.venv\Lib\site-packages\sklearn\metrics\pairwise.py:2466: DataConversionWarning: Data was converted to boolean for metric jaccard
  warnings.warn(msg, DataConversionWarning)
C:\Users\stas69\Documents\Studia\podstawy-uczenia-maszynowego-24-25\.venv\Lib\site-packages\sklearn\metrics\pairwise.py:2466: DataConversionWarning: Data was converted to boolean for metric jaccard
  warnings.warn(msg, DataConversionWarning)

Scikit-learn training time: 0.00
PyNNDescent training time: 5.43

Scikit-learn prediction time: 68.09
PyNNDescent prediction time: 0.64

Scikit-learn AUROC: 57.21%
PyNNDescent AUROC: 56.81%

Uważam że warto wykonać taką aproksymację poniważ mimo minimalnie gorszych wyników oraz dłuższego czasu treningu czas predykcji znacznie się skraca.

Wyszukiwanie z pomocą najbliższych sąsiadów

Naszym ostatnim zastosowaniem metod opartych o sąsiedztwo będzie wyszukiwanie przedmiotów. Konkretnie stworzymy wyszukiwarkę obrazów na podstawie wpisanego hasła.

Jest to zadanie multimodalne (multimodal), bo operujemy na dwóch różnych modalnościach: obrazie i tekście. Takie wyszukiwarki działają następująco:

1. W fazie treningu przedmioty (obrazy) są wektoryzowane siecią neuronową i przechowywane w indeksie.
2. Użytkownik wpisuje tekst, który jest wektoryzowany za pomocą tej samej sieci neuronowej, dlatego musi być ona multimodalna.
3. Wyszukujemy najbliższych sąsiadów (obrazy-wektory) dla tekstu-wektora.

Najbardziej znanym modelem multimodalnym jest CLIP, stworzony przez OpenAI. Jego architektura i wagi są dostępne publicznie, w przeciwieństwie do nowych rozwiązań tej firmy. Jest to model zasadniczo niewymagający fine-tuningu, działa bardzo dobrze w podstawowej formie dla wielu zadań. Został też wykorzystany jako element modelu DALL-E 2.

Nie musisz szczegółowo wiedzieć, jak działa ten model, ale pewne podstawy są bardzo ciekawe i przydatne. Wysokopoziomowo architekturę dobrze opisuje ten artykuł. Dla zainteresowanych cały artykuł, znacznie bardziej szczegółowy. CLIP składa się z dwóch enkoderów: dla tekstu i dla obrazu. Są to zwykłe pre-trenowane sieci, a konkretnie:

1. Dla obrazu - Visual Transformer (ViT). Jest to architektura, która używa transformerów i atencji do przetwarzania obrazów. Ideą jest, żeby pociąć obraz na sąsiadujące kawałki (patches), i traktować je jak słowa dla atencji.
2. Dla tekstu - Transformer, z pewnymi modyfikacjami.

Każda z tych sieci wektoryzuje wejście, dając embedding.

Zbiorem danych była kolekcja 400 milionów par (obraz, tekst), gdzie tekst opisywał zawartość obrazu. Następnie wykonano constrastive pre-training - model dostaje batch $N$ par (obraz, tekst) i ma przewidzieć, które z $N \times N$ możliwych kombinacji (obraz, tekst) faktycznie wystąpiły. W tym celu maksymalizuje podobieństwo cosinusowe par, które wystąpiły w batchu, a minimalizuje podobieństwo dla $N^2 - N$ możliwych kombinacji, których nie było w batchu. W tym celu jest używana entropia krzyżowa z pewnymi modyfikacjami.

Widać więc, że model jest nienadzorowany, a zbiór danych był ogromny, więc fine-tuning nie jest nam zasadniczo potrzebny. Można go użyć do klasyfikacji, ale tym nie będziemy się zajmować. Dodatkowo odpowiednią metryką odległości dla najbliższych sąsiadów będzie podobieństwo cosinusowe.

Biblioteka Sentence-Transformers powstała dla różnych modeli transformerowych do zadań nienadzorowanych, ale dodano do niej także bardzo wygodny interfejs dla modelu CLIP - przykład użycia.

Jako zbiór wykorzystamy Amazon Berkeley Objects (ABO) Dataset - zbiór obrazów produktów Amazona, stworzony we współpracy z University of California, Berkeley. Wykorzystamy zminiaturyzowany zbiór, w którym obrazy mają rozmiar 256x256 pikseli, bo na nasze potrzeby jest zupełnie wystarczający.

Zadanie 6 (3 punkty)

Wyszukiwarka z ANN

1. Ściągnij zbiór danych. Rozpakuj go obok tego notebooka. W katalogu images/metadata znajdziesz plik images.csv.gz - zawiera on ścieżki do obrazów.
2. Uzupełnij kod funkcji vectorize_images, która oblicza embeddingi obrazów za pomocą modelu CLIP. Przyda ci się ten tutorial. Ustaw wartość MAX_IMAGES tak, żeby wystarczyło ci RAMu i żeby proces trwał rozsądną ilość czasu (ale co najmniej kilka-kilkanaście minut). start_idx oraz end_idx to indeksy wierszy, na których w danej chwili operujemy - po prostu przesuwamy się o BATCH_SIZE w pętli for.
3. Uzupełnij kod klasy ImageSearch:

• w konstruktorze stwórz indeks za pomocą klasy NNDescent, korzystając z metryki cosinusowej i random_state=0, może się też przydać n_jobs=-1
• w metodzie query zwektoryzuj zapytanie tekstowe, a następ wyszukaj indeksy najbliższych sąsiadów-obrazów w indeksie (zwróć uwagę na to, co zwraca metoda .query() dla NNDescent), a następnie zwróć ścieżki do ich obrazów; zwróć uwagę na to, że .query() wymaga 2-wymiarowego wejścia, więc w naszym wypadku (1, 512)
• w metodzie show_images załaduj i wyświetl obrazy; może się przydać ten przykład

4. Przetestuj wyszukiwarkę kilkoma hasłami (metoda .search()). Prawdopodobnie najlepiej zadziałają hasła odpowiadające typowym produktom sprzedawanym w Amazonie.
5. Jakie widzisz zalety i wady takiego podejścia do wyszukiwania?

Uwaga: jeżeli masz mało RAMu lub pracujesz na Google Colab, to warto zrestartować Jupyter Notebooka, zamykając kernel. Zwolni ci to pamięć. Możesz też usunąć macierz z embeddingami (del embeddings) kiedy stworzysz indeks w konstruktorze ImageSearch.

!wget https://amazon-berkeley-objects.s3.amazonaws.com/archives/abo-images-small.tar

--2025-04-10 10:08:54--  https://amazon-berkeley-objects.s3.amazonaws.com/archives/abo-images-small.tar
Resolving amazon-berkeley-objects.s3.amazonaws.com (amazon-berkeley-objects.s3.amazonaws.com)... 3.5.21.80, 52.217.48.60, 54.231.139.17, ...
Connecting to amazon-berkeley-objects.s3.amazonaws.com (amazon-berkeley-objects.s3.amazonaws.com)|3.5.21.80|:443... connected.
HTTP request sent, awaiting response... 200 OK
Length: 3253381120 (3.0G) [application/x-tar]
Saving to: ‘abo-images-small.tar’

abo-images-small.ta 100%[===================>]   3.03G  36.6MB/s    in 1m 41s  

2025-04-10 10:10:36 (30.7 MB/s) - ‘abo-images-small.tar’ saved [3253381120/3253381120]

!tar -xf abo-images-small.tar

# prepare pd.Series with paths to images
image_paths = pd.read_csv("images/metadata/images.csv.gz")

# "height" and "width" are original image size, let's select
# only large ones - they are probably more interesting
mask = (image_paths["height"] >= 1000) & (image_paths["width"] >= 1000)
image_paths = image_paths.loc[mask, :]

# remove columns, leaving only pd.Series with paths
image_paths = image_paths["path"]
image_paths = "images/small/" + image_paths.astype(str)

image_paths

	path
554	images/small/46/4689906d.png
6614	images/small/0c/0cd7596c.jpg
6617	images/small/e6/e602a9af.jpg
6621	images/small/e7/e7cfeb02.jpg
6627	images/small/68/6839db4e.jpg
...	...
398207	images/small/6d/6d49d130.jpg
398208	images/small/b1/b163e0ea.jpg
398209	images/small/a1/a116d9d1.jpg
398210	images/small/9c/9c3e1158.jpg
398211	images/small/cf/cf112e38.jpg

342878 rows × 1 columns

dtype: object

from itertools import islice
from PIL import Image
import joblib
import torch
from sentence_transformers import SentenceTransformer
from tqdm.notebook import tqdm


MAX_IMAGES = 100_000
BATCH_SIZE = joblib.cpu_count(only_physical_cores=True)


def vectorize_images(image_paths: pd.Series) -> np.ndarray:
    device = "cuda" if torch.cuda.is_available() else "cpu"
    model = SentenceTransformer("clip-ViT-B-32", device=device)

    # batch iterator, based on Python itertools example
    def batched(iterable, n: int):
        it = iter(iterable)
        while batch := tuple(islice(it, n)):
            yield batch

    # CLIP embeddings have 512 dimensions
    embeddings = np.empty((MAX_IMAGES, 512))

    # iterate with tqdm, it will give us a nice progress bar
    with tqdm(total=MAX_IMAGES) as pbar:
        start_idx = 0
        for batch in batched(image_paths.iloc[:MAX_IMAGES], BATCH_SIZE):
            # update end_idx; remember to stop at MAX_IMAGES!
            end_idx = min(MAX_IMAGES, start_idx + BATCH_SIZE)
            # load images with Image.open()
            images = [Image.open(path) for path in batch]
            # calculate embeddings
            img_embeded = model.encode(images, convert_to_numpy=True)
            # save embeddings in the embeddings array
            embeddings[start_idx:end_idx] = img_embeded
            # update start_idx
            start_idx = end_idx
            pbar.update(BATCH_SIZE)

    return embeddings

embeddings = vectorize_images(image_paths)

/usr/local/lib/python3.11/dist-packages/huggingface_hub/utils/_auth.py:94: UserWarning: 
The secret `HF_TOKEN` does not exist in your Colab secrets.
To authenticate with the Hugging Face Hub, create a token in your settings tab (https://huggingface.co/settings/tokens), set it as secret in your Google Colab and restart your session.
You will be able to reuse this secret in all of your notebooks.
Please note that authentication is recommended but still optional to access public models or datasets.
  warnings.warn(

modules.json:   0%|          | 0.00/122 [00:00<?, ?B/s]

config_sentence_transformers.json:   0%|          | 0.00/116 [00:00<?, ?B/s]

README.md:   0%|          | 0.00/1.91k [00:00<?, ?B/s]

Xet Storage is enabled for this repo, but the 'hf_xet' package is not installed. Falling back to regular HTTP download. For better performance, install the package with: `pip install huggingface_hub[hf_xet]` or `pip install hf_xet`
WARNING:huggingface_hub.file_download:Xet Storage is enabled for this repo, but the 'hf_xet' package is not installed. Falling back to regular HTTP download. For better performance, install the package with: `pip install huggingface_hub[hf_xet]` or `pip install hf_xet`

merges.txt:   0%|          | 0.00/525k [00:00<?, ?B/s]

vocab.json:   0%|          | 0.00/961k [00:00<?, ?B/s]

preprocessor_config.json:   0%|          | 0.00/316 [00:00<?, ?B/s]

tokenizer_config.json:   0%|          | 0.00/604 [00:00<?, ?B/s]

special_tokens_map.json:   0%|          | 0.00/389 [00:00<?, ?B/s]

pytorch_model.bin:   0%|          | 0.00/605M [00:00<?, ?B/s]

config.json:   0%|          | 0.00/4.03k [00:00<?, ?B/s]

Using a slow image processor as `use_fast` is unset and a slow processor was saved with this model. `use_fast=True` will be the default behavior in v4.52, even if the model was saved with a slow processor. This will result in minor differences in outputs. You'll still be able to use a slow processor with `use_fast=False`.

  0%|          | 0/100000 [00:00<?, ?it/s]

from pynndescent import PyNNDescentTransformer, NNDescent
from mpl_toolkits.axes_grid1 import ImageGrid


class ImageSearch:
    def __init__(self, image_paths: pd.Series, embeddings: np.ndarray):
        # change to Numpy array to avoid .iloc[] and just index with []
        self.image_paths = image_paths.values

        device = "cuda" if torch.cuda.is_available() else "cpu"
        self.model = SentenceTransformer("clip-ViT-B-32", device=device)

        # create PyNNDescent index
        self.index = NNDescent(embeddings, metric="cosine", n_jobs=-1, random_state=0)

    def search(self, text: str, n_neighbors: int = 9) -> None:
        image_paths = self.query(text, n_neighbors)
        self.show(image_paths)

    def query(self, text: str, n_neighbors: int = 9) -> list[str]:
        text_emb = self.model.encode(text)
        text_emb = text_emb[np.newaxis, :]
        matching_indices = self.index.query(text_emb, k=n_neighbors)[0]
        return [self.image_paths[i] for i in matching_indices][0]

    def show(self, image_paths: list[str]) -> None:
        fig = plt.figure(figsize=(15.0, 15.0))
        grid = ImageGrid(
            fig,
            111,
            nrows_ncols=(3, 3),
            axes_pad=0.1,
        )

        for ax, im in zip(grid, image_paths):
            ax.imshow(Image.open(im))

        plt.show()

image_search = ImageSearch(image_paths, embeddings)

image_search.search("blue pillow")

image_search.search("PS4")

image_search.search("brown toy bear")

image_search.search("beer")

Zalety:

• Szybkie wyszukiwanie obrazów
• W większości wynikowe obrazy są trafne

Wady:

• przy niektórych opisach model skupia się osobno na poszczególnych elementach zapytania i zamiast zwrócić brązowego zabawkowego misia zwraca produkty które po prostu są brązowe
• W zapytaniu o PS4 model zwraca też inne konsole (switch)
• Przetworzenie tak dużej ilości obrazów jest kosztowne obliczeniowo i czasowo.

Zadanie dodatkowe (3 punkty)

Zapoznaj się z algorytmem Condensed Nearest Neighbors . Jest to algorytm typu training set reduction, który zmniejsza rozmiar zbioru treningowego do klasyfikacji, usuwając z niego "trywialne" próbki. Idea jest tutaj taka, że znaczenie mają zasadniczo tylko przykłady blisko granicy decyzyjnej, a te daleko od niej są oczywiste i można by tam ograniczyć liczbę próbek, przyspieszając algorytmy najbliższych sąsiadów.

Zaimplementuj algorytm Condensed Nearest Neighbors w postaci klasy kompatybilnej ze Scikit-learn, z metodami .fit() oraz .transform(). Wejściem powinny być macierz X oraz wektor y, a wyjściem zredukowane macierze, ograniczone tylko do interesujących próbek. Niestety ten algorytm jest czuły na kolejność przetwarzania próbek, dlatego ma element losowy - zaimplementuj argument random_state w konstruktorze i go użyj do przemieszania losowo zbioru (mogą się przydać generatory losowe z Numpy'a). Dodatkowo konstruktor powinien przyjmować argument metric (analogicznie do KNeighborsClassifier) i wykorzystywać daną metrykę jako sposób pomiaru odległości.

Przetestuj algorytm na obu zbiorach używanych w tym laboratorium:

1. Sprawdź, o ile próbek udało się zredukować zbiory treningowe.
2. Porównaj szybkość treningu, szybkość predykcji oraz jakość wyników. Zmierz też czas wykonywania samego algorytmu Condensed Nearest Neighbors.
3. Skomentuj, czy twoim zdaniem warto korzystać z takiego rozwiązania.

Przydatne źródła:

• dykusja na DataScience Stack Exchange
• prezentacja University of Leicester o kNN i CNN

Dla zainteresowanych polecam artykuły naukowe:

Hart, Peter. "The condensed nearest neighbor rule (corresp.)." IEEE transactions on information theory 14.3 (1968): 515-516. link

Angiulli, Fabrizio. "Fast condensed nearest neighbor rule." Proceedings of the 22nd international conference on Machine learning. 2005. link